The binomial coefficients nCk can be arranged in triangular form, Pascal's triangle, like this:
1
1 1
1 2 1
1 3 3 1
1 4 6 4 1
1 5 10 10 5 1
1 6 15 20 15 6 1
1 7 21 35 35 21 7 1
.........
It can be seen that the first eight rows of Pascal's triangle contain twelve distinct numbers:
1, 2, 3, 4, 5, 6, 7, 10, 15, 20, 21 and 35.
A positive integer n is called squarefree if no square of a prime divides n. Of the twelve distinct
numbers in the first eight rows of Pascal's triangle, all except 4 and 20 are squarefree. The sum of
the distinct squarefree numbers in the first eight rows is 105.
Find the sum of the distinct squarefree numbers in the first 51 rows of Pascal's triangle.
In [1]:
N = 51
def binom(n, k):
if k > n // 2:
k = n - k
b = 1
for i in range(k):
b = (b * (n - i)) // (i + 1)
return b
In [2]:
from sympy import primerange, nextprime
MAXIMUM_ROW_NUMBER = 100
smallprimes = list(primerange(2, MAXIMUM_ROW_NUMBER + 1))
In [3]:
def squarefree_binom(n, k):
assert n <= MAXIMUM_ROW_NUMBER
for p in smallprimes:
if p > n:
return True
if binom_order(n, k, p) > 1:
return False
return True
The exponent of a prime $p$ in the prime factorization of $n!$ is equal to $(n-s)/(p-1)$, where $s$ is the sum of the base-$p$ digits of $n$. This is a consequence of Legendre's formula. We can use this formula to calculate the exponent of $p$ in the prime factorization of a binomial coefficient.
In [4]:
def sum_of_digits(n, p):
return n if n < p else (n % p) + sum_of_digits(n // p, p)
def fact_order(n, p):
return (n - sum_of_digits(n, p)) // (p - 1)
def binom_order(n, k, p):
return fact_order(n, p) - fact_order(k, p) - fact_order(n-k, p)
In [5]:
results = set()
for n in range(N):
for k in range(n//2 + 1):
if squarefree_binom(n, k):
results.add(binom(n, k))
print(sum(results))